##### #############################################################
#####
#####   Input: Clean experimental data
#####   Output: Analysis and figures for the main paper
#####

setwd("/Users/lotte/Dropbox/projects/productivity_experiment/replication")

rm(list = ls())

library(data.table) # CRAN v1.14.2
library(plyr) # CRAN v1.8.6
library(dplyr) # CRAN v1.0.9
library(tidyverse) # CRAN v1.3.1
library(cjoint) # CRAN v2.1.0
library(ggplot2) # CRAN v3.3.6
library(broom) # CRAN v0.7.11
library(patchwork) # CRAN v1.1.1
library(estimatr) # CRAN v0.30.6

# Load data

load("data/productivity.Rdata")

# Analysis

results_amce <-
  lm_robust(
    preference ~ mp_gender + mp_party + committee_membership +
      issue_campaigning + voting_legislation + constituency_responsiveness,
    data = productivity,
    clusters = ID,
    se_type = "stata"
  )

amce <- tidy(results_amce, conf.int = T)[2:7, ]

amce$attributes <-
  factor(
    NA,
    levels = c(
      "MP Gender",
      "MP Party",
      "Committee Membership",
      "Constituency Responsive",
      "Campaign Successful",
      "Voting Productive"
    )
  )
amce$attributes[amce$term == "mp_genderWoman"] <- "MP Gender"
amce$attributes[amce$term == "mp_partyLabour"] <- "MP Party"
amce$attributes[amce$term == "committee_membershipSits on several"] <-
  "Committee Membership"
amce$attributes[amce$term == "constituency_responsivenessOften"] <-
  "Constituency Responsive"
amce$attributes[amce$term == "issue_campaigningSuccessfully campaigned"] <-
  "Campaign Successful"
amce$attributes[amce$term == "voting_legislationMore productive"] <-
  "Voting Productive"
table(amce$attributes)

figure_2_amce <-
  ggplot(amce, aes(x = estimate, y = attributes)) + geom_point(size = 3) + ylab("") + xlab("") +
  geom_linerange(aes(xmin = conf.low, xmax = conf.high)) +
  theme_bw() +
  geom_vline(xintercept = 0, linetype = 2) +
  xlim(-0.1, 0.3) +
  theme(axis.text = element_text(size = 11),
        axis.title = element_text(size = 11)) +
  scale_y_discrete(
    limits = rev,
    labels = c(
      "Voting\n Productive",
      "Campaign\n Succcessful",
      "Constituency\n Responsive",
      "Committee\n Member",
      "Conservative\n MP",
      "Woman\n MP"
    )
  ) +
  labs(x = "Change in Pre(Preferred Member of Parliament)", y = "", title = "")
ggsave("plots/figure_2_amce.png",
       figure_2_amce,
       height = 4,
       width = 6)

performance_bivariate <-
  lm_robust(
    perceived_performance ~ mp_gender + mp_party + committee_membership +
      issue_campaigning + voting_legislation + constituency_responsiveness,
    data = productivity,
    clusters = ID,
    se_type = "stata"
  )

electability_bivariate <-
  lm_robust(
    electability  ~ mp_gender + mp_party + committee_membership +
      issue_campaigning + voting_legislation + constituency_responsiveness,
    data = productivity,
    clusters = ID,
    se_type = "stata"
  )

unconditional <- bind_rows(
  tidy(performance_bivariate, conf.int = T)[2:7, ] %>% mutate(outcome = "Job Performance"),
  tidy(electability_bivariate, conf.int = T)[2:7, ] %>% mutate(outcome =
                                                                 "Electability")
)

unconditional$attributes <-
  factor(
    NA,
    levels = c(
      "MP Gender",
      "MP Party",
      "Committee Membership",
      "Constituency Responsive",
      "Campaign Successful",
      "Voting Productive"
    )
  )
unconditional$attributes[unconditional$term == "mp_genderWoman"] <-
  "MP Gender"
unconditional$attributes[unconditional$term == "mp_partyLabour"] <-
  "MP Party"
unconditional$attributes[unconditional$term == "committee_membershipSits on several"] <-
  "Committee Membership"
unconditional$attributes[unconditional$term == "constituency_responsivenessOften"] <-
  "Constituency Responsive"
unconditional$attributes[unconditional$term == "issue_campaigningSuccessfully campaigned"] <-
  "Campaign Successful"
unconditional$attributes[unconditional$term == "voting_legislationMore productive"] <-
  "Voting Productive"
table(unconditional$attributes)

figure_3 <-
  ggplot(unconditional, aes(x = estimate, y = attributes)) + geom_point(size = 3) + ylab("") + xlab("") +
  geom_linerange(aes(xmin = conf.low, xmax = conf.high)) +
  theme_bw() +
  facet_wrap( ~ outcome, scale = "free") +
  geom_vline(xintercept = 0, linetype = 2) +
  theme(axis.text = element_text(size = 11),
        strip.text.x = element_text(size = 12)) +
  xlim(-0.5, 1.6) +
  scale_y_discrete(
    limits = rev,
    labels = c(
      "Voting\n Productive",
      "Campaign\n Succcessful",
      "Constituency\n Responsive",
      "Committee\n Member",
      "Conservative\n MP",
      "Woman\n MP"
    )
  )
ggsave(
  "plots/figure_3_unconditional_evaluations.png",
  figure_3,
  width = 7,
  height = 5
)

# Conditional by MP gender (binary productivity measures)

electability_levels_woman <-
  lm_robust(
    electability ~ committee_membership + issue_campaigning + voting_legislation +
      constituency_responsiveness,
    data = productivity,
    subset = mp_gender == "Woman",
    clusters = ID,
    se_type = "stata"
  )
electability_levels_man <-
  lm_robust(
    electability ~ committee_membership + issue_campaigning + voting_legislation +
      constituency_responsiveness,
    data = productivity,
    subset = mp_gender == "Man",
    clusters = ID,
    se_type = "stata"
  )
performance_levels_woman <-
  lm_robust(
    perceived_performance ~ committee_membership + issue_campaigning + voting_legislation +
      constituency_responsiveness,
    data = productivity,
    subset = mp_gender == "Woman",
    clusters = ID,
    se_type = "stata"
  )
performance_levels_man <-
  lm_robust(
    perceived_performance ~ committee_membership + issue_campaigning + voting_legislation +
      constituency_responsiveness,
    data = productivity,
    subset = mp_gender == "Man",
    clusters = ID,
    se_type = "stata"
  )
preference_levels_woman <-
  lm_robust(
    preference ~ committee_membership + issue_campaigning + voting_legislation +
      constituency_responsiveness,
    data = productivity,
    subset = mp_gender == "Woman",
    clusters = ID,
    se_type = "stata"
  )
preference_levels_man <-
  lm_robust(
    preference ~ committee_membership + issue_campaigning + voting_legislation +
      constituency_responsiveness,
    data = productivity,
    subset = mp_gender == "Man",
    clusters = ID,
    se_type = "stata"
  )

attributes_conditional <-
  bind_rows(
    tidy(electability_levels_woman, conf.int = T)[2:5, ] %>% mutate(gender =
                                                                      "Woman", outcome = "Electability"),
    tidy(electability_levels_man, conf.int =
           T)[2:5, ] %>% mutate(gender = "Man", outcome = "Electability"),
    tidy(performance_levels_woman, conf.int = T)[2:5, ] %>%
      mutate(gender = "Woman", outcome = "Job Performance"),
    tidy(performance_levels_man, conf.int = T)[2:5, ] %>%
      mutate(gender = "Man", outcome = "Job Performance"),
    tidy(preference_levels_woman, conf.int = T)[2:5, ] %>%
      mutate(gender = "Woman", outcome = "MP Preference"),
    tidy(preference_levels_man, conf.int = T)[2:5, ] %>%
      mutate(gender = "Man", outcome = "MP Preference")
  )

figure_4 <-
  ggplot(attributes_conditional,
         aes(
           x = estimate,
           y = term,
           group = gender,
           color = gender
         )) +
  geom_point(position = position_dodge(width = 0.2), size = 2) + ylab("") + xlab("") +
  theme_bw() +
  geom_linerange(position = position_dodge(width = 0.2),
                 aes(
                   xmin = conf.low,
                   xmax = conf.high,
                   linetype = gender
                 )) +
  facet_wrap( ~ outcome, scales = "free") +
  scale_color_manual(values = c("lightseagreen", "grey40")) +
  labs(color = "MP gender",
       linetype = "MP gender",
       shape = "MP gender") +
  geom_vline(xintercept = 0, lty = "dashed") +
  theme(
    axis.text = element_text(size = 9),
    strip.text.y = element_text(size = 9),
    strip.text.x = element_text(size = 10),
    legend.title = element_text(size = 9)
  ) +
  scale_y_discrete(
    labels = c(
      "Committee \n Member",
      "Constituency \n Responsive",
      "Campaign \n Succcessful",
      "Voting \n Productive"
    )
  )
ggsave(
  "plots/figure_4_conditional_binary.png",
  figure_4,
  width = 10,
  height = 4
)

# Conditional by MP gender (continuous productivity measure)

model_1 <- lm_robust(
  preference ~ objective_performance * mp_gender,
  data = productivity,
  clusters = ID,
  se_type = "stata"
)
model_2 <-
  lm_robust(
    electability ~ objective_performance * mp_gender,
    data = productivity,
    clusters = ID,
    se_type = "stata"
  )
model_3 <-
  lm_robust(
    perceived_performance ~ objective_performance * mp_gender,
    data = productivity,
    clusters = ID,
    se_type = "stata"
  )

new_data_man <- data.frame(mp_gender = "Man",
                           objective_performance = c(0:4))
new_data_woman <- data.frame(mp_gender = "Woman",
                             objective_performance = c(0:4))

model_1_man <-
  predict(model_1, interval = "confidence", newdata = new_data_man)
model_1_woman <-
  predict(model_1, interval = "confidence", newdata = new_data_woman)
model_2_man <-
  predict(model_2, interval = "confidence", newdata = new_data_man)
model_2_woman <-
  predict(model_2, interval = "confidence", newdata = new_data_woman)
model_3_man <-
  predict(model_3, interval = "confidence", newdata = new_data_man)
model_3_woman <-
  predict(model_3, interval = "confidence", newdata = new_data_woman)

model_1_for_plot <-
  as.data.frame(rbind(model_1_man$fit, model_1_woman$fit)) %>%
  mutate(
    gender = c(
      "Man",
      "Man",
      "Man",
      "Man",
      "Man",
      "Woman",
      "Woman",
      "Woman",
      "Woman",
      "Woman"
    ),
    outcome = c("MP preference"),
    objective_performance = c(0, 1, 2, 3, 4, 0, 1, 2, 3, 4)
  )

model_2_for_plot <-
  as.data.frame(rbind(model_2_man$fit, model_2_woman$fit)) %>%
  mutate(
    gender = c(
      "Man",
      "Man",
      "Man",
      "Man",
      "Man",
      "Woman",
      "Woman",
      "Woman",
      "Woman",
      "Woman"
    ),
    outcome = c("Electability"),
    objective_performance = c(0, 1, 2, 3, 4, 0, 1, 2, 3, 4)
  )

model_3_for_plot <-
  as.data.frame(rbind(model_3_man$fit, model_3_woman$fit)) %>%
  mutate(
    gender = c(
      "Man",
      "Man",
      "Man",
      "Man",
      "Man",
      "Woman",
      "Woman",
      "Woman",
      "Woman",
      "Woman"
    ),
    outcome = c("Job Performance"),
    objective_performance = c(0, 1, 2, 3, 4, 0, 1, 2, 3, 4)
  )

models_for_plot <- rbind(model_1_for_plot,
                         model_2_for_plot,
                         model_3_for_plot)

figure_5 <-
  ggplot(models_for_plot,
         aes(y = fit, x = objective_performance, group = gender)) +
  geom_ribbon(aes(
    ymin = lwr,
    ymax = upr,
    linetype = gender,
    fill = gender
  ), alpha = 0.3) +
  theme_bw() +
  scale_fill_manual(values = c("lightseagreen", "grey40")) +
  scale_colour_manual(values = c("lightseagreen", "grey40")) +
  geom_line(aes(
    group = gender,
    color = gender,
    linetype = gender
  )) +
  labs(color = "MP gender",
       linetype = "MP gender",
       fill = "MP gender") +
  facet_wrap( ~ outcome, scales = "free") +
  ylab("") + xlab("Objective performance") +
  theme(
    axis.text = element_text(size = 11),
    strip.text.x = element_text(size = 10),
    legend.title = element_text(size = 10),
    axis.title = element_text(size = 11)
  )
ggsave(
  "plots/figure_5_conditional_continuous.png",
  figure_5,
  width = 10,
  height = 5
)

## Party congruence -- continuous

model_6 <-
  lm_robust(
    preference ~ objective_performance * party_congruence,
    data = productivity,
    clusters = ID,
    se_type = "stata"
  )
model_7 <-
  lm_robust(
    perceived_performance ~ objective_performance * party_congruence,
    data = productivity,
    clusters = ID,
    se_type = "stata"
  )
model_8 <-
  lm_robust(
    electability ~ objective_performance * party_congruence,
    data = productivity,
    clusters = ID,
    se_type = "stata"
  )

new_data_congruent <- data.frame(party_congruence = "TRUE",
                                 objective_performance = c(0:4))
new_data_incongruent <- data.frame(party_congruence = "FALSE",
                                   objective_performance = c(0:4))
model_6_congruent <-
  predict(model_6, interval = "confidence", newdata = new_data_congruent)
model_6_incongruent <-
  predict(model_6, interval = "confidence", newdata = new_data_incongruent)

model_6_congruent_for_plot <-
  as.data.frame(rbind(model_6_congruent$fit, model_6_incongruent$fit)) %>%
  mutate(
    congruence = c(
      "Congruent",
      "Congruent",
      "Congruent",
      "Congruent",
      "Congruent",
      "Incongruent",
      "Incongruent",
      "Incongruent",
      "Incongruent",
      "Incongruent"
    ),
    outcome = c("MP Preference"),
    objective_performance = c(0, 1, 2, 3, 4, 0, 1, 2, 3, 4)
  )

model_7_congruent <-
  predict(model_7, interval = "confidence", newdata = new_data_congruent)
model_7_incongruent <-
  predict(model_7, interval = "confidence", newdata = new_data_incongruent)

model_7_congruent_for_plot <-
  as.data.frame(rbind(model_7_congruent$fit, model_7_incongruent$fit)) %>%
  mutate(
    congruence = c(
      "Congruent",
      "Congruent",
      "Congruent",
      "Congruent",
      "Congruent",
      "Incongruent",
      "Incongruent",
      "Incongruent",
      "Incongruent",
      "Incongruent"
    ),
    outcome = c("Job Performance"),
    objective_performance = c(0, 1, 2, 3, 4, 0, 1, 2, 3, 4)
  )

model_8_congruent <-
  predict(model_8, interval = "confidence", newdata = new_data_congruent)
model_8_incongruent <-
  predict(model_8, interval = "confidence", newdata = new_data_incongruent)

model_8_congruent_for_plot <-
  as.data.frame(rbind(model_8_congruent$fit, model_8_incongruent$fit)) %>%
  mutate(
    congruence = c(
      "Congruent",
      "Congruent",
      "Congruent",
      "Congruent",
      "Congruent",
      "Incongruent",
      "Incongruent",
      "Incongruent",
      "Incongruent",
      "Incongruent"
    ),
    outcome = c("Electability"),
    objective_performance = c(0, 1, 2, 3, 4, 0, 1, 2, 3, 4)
  )

party_congruence <-
  rbind(
    model_6_congruent_for_plot ,
    model_7_congruent_for_plot,
    model_8_congruent_for_plot
  )

figure_6 <-
  ggplot(party_congruence,
         aes(y = fit, x = objective_performance, group = congruence)) +
  geom_ribbon(aes(
    ymin = lwr,
    ymax = upr,
    linetype = congruence,
    fill = congruence
  ),
  alpha = 0.3) +
  theme_bw() +
  scale_fill_manual(values = c("lightseagreen", "grey40")) +
  scale_colour_manual(values = c("lightseagreen", "grey40")) +
  geom_line(aes(
    group = congruence,
    color = congruence,
    linetype = congruence
  )) +
  labs(color = "Party congruence",
       linetype = "Party congruence",
       fill = "Party congruence") +
  facet_wrap( ~ outcome, scales = "free") +
  ylab("") + xlab("Objective performance") +
  theme(
    axis.text = element_text(size = 11),
    strip.text.x = element_text(size = 10),
    legend.title = element_text(size = 10),
    axis.title = element_text(size = 11)
  )
ggsave(
  "plots/figure_6_party_congruence.png",
  figure_6,
  width = 10,
  height = 5
)

